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ABSTRACT 



A numerical model for heat transfer in laminar duct flows has been developed 
using the finite difference method to explore the significance and extent of “back- 
conduction” at low Peclet numbers. The calculations have been carried out for flows 
between parallel plates and in circular tubes by using different Peclet numbers in the 
range of 0.05 to 100. For both situations constant heat flux and constant wall temperature 
boundary conditions were used. The validity of the results has been checked by 
comparison with some existing results in the literature, and extended to a wider range of 
parameters including conjugate wall conduction effects. The results are presented for bulk 
mean temperature variation, Nusselt number behavior, and energy absorbed before the 
heated section, for cases with and without wall conduction. Such axial conduction effects 
may be an important feature in the thermal characterization of microtubes, which are to 
be used in microheat exchangers. 
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LIST OF SYMBOLS 



SYMBOL 


DESCRIPTION 


UNITS 


Cp 


Specific heat at constant pressure 


J/kg.K 


d 


The diameter of the circular tube 


m 


Dh 


Hydraulic diameter for parallel plates, Dh=4h 


m 


^absorbed 


Absorbed energy inside the fluid from the entrance to x=0 


Joule 


E* 


Dimensionless absorbed energy 


unity 


h 


Convection heat transfer coefficient 


W/mlK 


h 


The distance from the centerline to the surface of the 






parallel plate 


m 


hx 


The distance between two grid points in x direction 


m 


hy 


The distance between two grid points in y direction 


m 


k 


Thermal conductivity 


W/m.K 


Nu 


Nusselt number 


unity 


Pe, PcDh 


Peclet number, Pe=Re.Pr 


unity 


Pr 


Prandtl number 


unity 


<io 


Heat flux from the surface in constant heat flux case 


wW 


r 


Distance from the center of the circular tube 


m 


ro 


Radius of the circular tube 


m 


+ 

r 


Dimensionless radius, rfro 


unity 


Re 


Reynolds number 


unity 


T,t 


Temperature 




^5 ^0 


Wall temperature 


^KX 


te 


Fluid entrance temperature 


^KX 




Mean temperature 


^KX 


U 


Velocity 


m/sec 


U 


Dimensionless velocity 


unity 


V 


Average velocity 


m/sec 


Vr 


Radial velocity of the fluid inside the circular tube 


m/sec 


W 


Relaxation parameter 


unity 


X 


Dimensionless distance in x direction 


unity 


y 


Distance from the centerline of the parallel plates 


m 


V 


Kinematic viscosity 


mVsec 


a 


Thermal diffusivity 


m^/sec 


M 


Viscosity 


kg/sec.m 


P 


Density 


kg/m^ 


6 


Bulk temperature 


unity 




Bulk mean temperature 


unity 
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I. INTRODUCTION 



The heat transfer problem of laminar fluid flow in ducts which is known as the 
Graetz problem has many applications in technology and has been studied extensively 
since Graetz (1885). The classical Graetz problem considers the forced convection heat 
transfer of the fluids flowing in ducts neglecting axial conduction effects in the fluid or 
the wall. 

This study is on the axial conduction effects of the flows whose Peclet numbers 
range is from those of liquid metals where axial conduction may not be neglected, to 
those of oils where the axial conduction has almost no effect on the temperature 
distribution inside the ducts. Developing micromachines nowadays and in the future, the 
axial conduction effects may be a very important feature in microheat exchangers. 

The effects of axial conduction on the Hagen-Poiseuille flows between parallel plates and 
in circular cylinders for constant wall temperature and constant heat flux have been 
studied in this thesis. The governing equations have been solved using the finite 
difference numerical scheme. The validity of the results has been checked by comparison 
with some existing results in the literature and extended to a wider range of parameters 
including conjugate wall conduction effects. The results are presented for bulk mean 
temperature variation, Nusselt number behavior, and energy absorbed before the heated 
section, for cases with and without wall conduction. 
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II. BACKGROUND 



A. HEAT TRANSFER IN LAMINAR DUCT FLOWS 

Consider a flow in a circular tube where fluid enters the tube with a uniform 
velocity. The viscous effects are important in this flow and a boundary layer develops as 
X increases. The flow is “fiilly developed” in the region where the velocity gradient 
doesn’t change anymore with increasing x (Figure 1). 









i r; Heal transfer by way of constant heat 

V' flux or constant wall temperature 




^ ryt 1''-^ 

ufrLuM ^ ^ t 




X 


LAMINAI 


1 FULLY DEVELOPED FLOW 

x = 0 



Figure 1. Laminar fully developed flow 

The Reynolds Number for the duct flows is defined as 

( 2 . 1 ) 

V 

where V is the mean fluid velocity, v is kinematic viscosity, v = —, and Dfj is the 

P 

hydraulic diameter, £, _ ^ Cross ^Sectional _Area . For the laminar duct flows, 

^ Wetted _Perimeter 

Re^ < 2300 . 

The Prandtl number relates the temperature distribution to the velocity distribution 
and defined as 

Pr = - (2.2) 

a 
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where a is the thermal diffusivity . If the Prandtl number is one, the velocity and the 
temperature profiles develop together and at the same rate. The Prandtl number is: 

Pr « 1 for the liquid metals 

Pr = 1 for the gases 

Pr»l for the oils 

For the liquid metals, where the Pr«l, the energy diffusion rate is much more 
than the momentum diffusion rate. It is opposite for the oils, in which Pr»l and the 
velocity profile develops faster than the temperature profile in this case. The value of the 
Prandtl number strongly affects the relative growth of the velocity and thermal boundary 
layers. 

Peclet number is defined as 

Pe=Re.Pr (2.3) 

and it is a measure of the quantity of the axial heat conduction effects in the fluid. The 
axial conduction is assumed negligible for Pe>10 and the axial conduction term can be 
assumed small in the governing equation for this case. But the axial conduction effects 
can be significant when the Peclet number is smaller. The purpose of this numerical study 
is to show the importance of the axial conduction effects and the heat absorption upto the 
entrance firom where the heating starts when the Peclet number is small. 

The Nusselt number is defined as 



Nu = ^ (2.4) 

k 

where k is the thermal conductivity of the fluid. It is a non-dimensional number where the 
convection coefficient h is calculated. 

The laminar duct flow being considered is assumed to be hydrodynamically fully 
developed before any heating effects are considered. The velocity profile in the laminar 
fully developed region is: 

For the circular tube: 



u(r) = 2V 



1 - 






v^oy 



(2.5) 
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For the parallel plates: 



u(y) = -v 







1 

H 






IhJJ 



( 2 . 6 ) 



where r and y are measured from the centerline, ro is tube radius and h is the distance 
from the centerline to the plates. [Ref.l] 

1. The Heat Transfer 

The heat transfer for a duct flow can be expressed using Newton’s law of cooling, 
q =h(Tw-Tm), where h is the convection heat transfer coefficient, Tw is the wall 
temperature and Tn, is the mean temperature of the fluid, where Tm is 



jpuc^TdA 

ip _ Area 



(2.7) 



me. 



For the constant Cv and incompressible flow through the circular tube, Tm is 

2 



T„ =■ 



Vr, 



— JuTrdr 

o 0 



( 2 . 8 ) 



where for the flow between the parallel plates is 



= iJuTdy 



(2.9) 



y and r are measured from the centerline, ro is the radius of the tube and h is the distance 
from the centerline to the plate. [Ref. 1] 

The energy equations for the cylindrical and rectangular coordinates are 



and 



di di 

ox or 



di di di 

up V\p V-wp 

dx dy dz 



’l d 


{kr^^ 


1 d 


1 




{k^\ 


r dr 


1 


^ r^ d(f> 


V^) 


1 dx 


1 ^)_ 



d (, dt 



dx 



\ 



f 



dt 



dy{ dyj 



+ — 
dz 



= 0 



f; 

k — 

I 



( 2 . 10 ) 



= 0 ( 2 . 11 ) 
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where ”i” is the enthalpy. These equations are applied to the problem including the 
boundary conditions and the temperature distribution is calculated analytically or 
numerically. [Ref. 2] 

2. The Thermally Fully Developed Flow 



In thermally folly developed case, the relative shape of the profile no longer 
changes although the temperature profile changes with x. The criterion for a folly 
developed temperature profile for a circular tube is 



T,„ -T^ 
dr 



= const 






( 2 . 12 ) 



This condition is reached in a duct flow whether the case is uniform heat flux (q 
through the wall is constant) or uniform wall temperature (The wall temperature is 
constant). These two cases have a lot of applications in engineering. For example, 
constant wall temperature in boiling or condensation, or constant heat flux using the 
electrical heater. In the thermally folly developed flow, the convection coefficient “h” is 
constant, independent of x. [Ref. 1] 



B. BACKGROUND STUDIES 



The heat transfer solution for laminar folly developed parabolic velocity profile 
flow inside a circular tube in the thermal region and subject to uniform tube wall 
temperature was treated for the first time by Graetz in 1883 and is known in the heat 
transfer literature as Graetz problem. In his problem Graetz neglected the axial 
conduction and solved the problem. [Ref. 3] 

The problem statement is; 



1 (. 2\d0, d^9, 1 dO. 

2 ' dx. dx: r. dr. 



(2.13) 



A = 0 
de ^ 

dr. 

9.=\ 



at r. = 1 
at r. = 0 

at X. = 1 
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where 






r. = • 



Tin -To 



X. - ■ 



Re^Pr 

and To is the wall temperature, Tin is the initial temperature, r is the radial distance, x is 
the axial distance, D is the diameter. Rep is the Reynolds number and Pr is the Prandtl 
number. [Ref. 4] 

The velocity profile is fully developed much before the temperature profile when 
the Prandtl number is high relative to one. For such situations, Graetz solution is 
reasonably well justified. But the Prandtl number is very small when the fluid is liquid 
metal and therefore the axial conduction effects may not be neglected. For such cases, 
since the axial conduction effects are not included to the result, the solution may give 
inaccurate results. 

Michelsen and Villadsen investigated the Graetz problem with axial heat 
conduction for circular tube constant wall temperature case by using a numerical 
procedure. They used a method that is the combination of orthogonal collocation and 
matrix diagonalization. They didn’t include the wall conduction effects in the problem. 

The problem may be defined as: 






X 



30 ^ 1 8^0 



+ ■ 



dy X dx\ dx) Pe^ dy 



(2.14) 



where 



k 



r z 

x = — , y = 

R PeR 



T-Tq ^ 









The boundary conditions are: 

y —> —00 0 = 1 

x = 0 ^ = 0 

dx 



x = 1 



^ = 0 
dx 

0 = 0 



for y < 0 
for y>0 
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where Vz is the average fluid velocity, R is the tube radius, z is the axial distance, r is the 
radial distance, Tb is the fluid temperature at z ^ -co , t is the bulk mean temperature 
and To is the wall temperature at z > 0 . They plotted the results, Nusselt number vs. y, 

heat flux vs. y and 0 vs. y where the axial heat conduction effects are easily seen. But 
they didn’t include the axial heat transfer inside the wall, which increases the fluid 
temperature very much upto the entrance as the wall conductivity gets larger. [Ref. 5] 

X. Yin and H. H. Bau included the axial conduction effects of the wall to the 
internal flow through circular tube, by using two parameters, duct’s outer/inner radius 
ratio and fluid/wall thermal conductivity ratio. They plotted the graphs of “temperature 
vs. radius of the circular tube” and “Nusselt number vs. Peclet number”, including the 
axial fluid and wall conduction effects. [Ref. 6] 

C. NUMERICAL METHOD 

The governing equations are steady state elliptic partial differential equations. The 
temperature T(x,y) throughout the domain must satisfy the governing equation and the 
boundary conditions along the entire boundary. The finite difference method was used in 
this numerical study to express the governing equations numerically. The central, forward 
and backward difference expressions used in this numerical study for the grid points in 
figure(2) are as follows : 




Figure (2). The grid points in finite difference methods 
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The central-difference expressions for the first and second derivatives with error 
order of h are: 



^ ki+lJ) ~ hi-lj) 
dx 2h^ 



(2.15a) 



dy 2hy 



(2.15b) 



dx^ {h^f 



(2.16a) 



_ ku+^) hu]-\) ^h'j) 

(hyf 

2 • 

The forward-difference expression with error order of h is: 

~ hiJ+2) + ^kiJ+d ~ 
dy~ 2hy 

2 

The backward-difference expression with error order of h is: 



(2.16b) 



(2.17) 



_ h‘J+2) ^h'J) 

^ 2hy 



[Ref. 7] 



(2.18) 



Substituting these finite-difference expressions into the governing equations, the problem 
is reduced to a set of linear algebraic equations, depending on the number of grid points 
used. Usually large number of grid points is desirable, but the number of equations to be 
solved becomes too large in this case. Gauss Seidel iteration method can be used to obtain 
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a solution. All the temperatures of the grid points are then make equal to the initial 
temperature at the beginning of the numerical solution. The temperatures at the grid 
points are calculated by using the governing equations and the boundary conditions. After 
the temperature of a grid point is calculated using the finite difference equations, this 
updated temperature is used for the next grid and the temperature values converge to the 
exact temperature a little more in every iteration. A technique called over relaxation is 
used to speed the convergence of the Gauss Seidel method when applied to the elliptic 
partial differential equaitions. This technique is: 



_new) ^ ( 2 ' 1 ^) 

w is known as the relaxation parameter. For overrelaxation w is between 1 and 2; for 
underrelaxation between 0 and 1. For this method, overrelaxation must be used to speed 
up the convergence. [Ref. 7] 

Forsythe and Wasow show that for a 45 by 45 mesh grid the optimum value of w 
is around 1.870. They also point out that when this relaxation parameter is used, the 
convergence is approximately 30 times faster than the usual process (w=l). [Ref 8] 
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III. GOVERNING EQUATIONS 
A. DERIVATION OF THE EQUATIONS 

1. Flow Inside the Circular Tube, Constant Wall Temperature 

Case 

For flow in a circular tube, the viscous energy equation in cylindrical coordinate 
system is 



di di 

up- + Vrp—- 

ox or 



1 



1 




1 a 

j 




a 

4_ 


(k^]] 


r dr ^ 


[ Sr] 


r^ d<p 


1 


dx 


[ dx)_ 



= 0 



(3.1) 



where di - cdt + — dP and is the radial velocity [Ref.l]. Then the equation becomes 
P 



for the constant k 



dt dt 

upc — + v,pc 

dx dr 



k d ( dt\ k dh , d^t 

+ -^ :r + k- 



r V dr 



= 0 



(3.2) 



df dx^ 

The temperature distribution is symmetric ip^tj d(f>^^=Q and the flow is 
hydrodynamically fully developed, v,. = 0 . Equation then becomes 



r dr \ dr 



+ - 



dx^ 



k Jdx 



(3.3) 



a = kj pc , the thermal diffiisivity of the fluid. So the governing equation for the circular 
duct flow is 



d^t 1 dt d^t u dt 



• + +■ 



dr^ r dr dx^ a dx 

Now, to non-dimensionalize the problem, assume that 



(3.4) 
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0=^ 



tn -t 



to-te 



■, r = — , u = — , X = 



+ _xIto 



RePr 



, =2{\-r* ) 



The governing equation becomes 









dx^ dr^ dr^ Pe^ dx^ 



(3.5) 



2. Flow Inside the Circular Tube, Constant Heat Flux Case 



The governing equation for the circular duct flow is (3.4) 



d^t 1 d~t _u dt 

dr^ r dr dx^ a dx 



Now, to non-dimensionalize the problem, assume that 



0= ^ 



qldjk 



u 

V 



+ 2x1 d 



Re^ Pr 



, u^=2i\-r^) 



The governing equation becomes 



, Ade d'^e 1 de i d'^e 

1-r^ P— = ^ + — - + ■ 



dx^ dr* r* dr* Pe^ dx* 



(3.6) 



3. Flow Between the Parallel Plates, Constant Wall Temperature 



Case 



The viscous energy equation in cartesian coordinate system for the flow between 
the parallel plates is 



di di di 

up hvp \-wp — 

dx dy dz 



dx 



(i ^0 

k — 



dy 



f, 5 /^ 
k — 



+ — 
dz 



(i ^0 

k — 

I 



= 0 (3.7) 
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where di = cdt + — dP and is the radial velocity [Ref. 1 ] . Then the equation becomes 
P 

for the constant k 



dt dt dt , 

dx ^ dy dz 



dx 



(dt\ d 
+ — 



\dxj 



dy 






d (dt^ 
+ — 



dz 



\dzj 



= 0 (3.8) 



3t 

The temperature distribution on z direction is constant, — = 0 and the flow is laminar and 

dz 

hydrodynamically fully developed, v = 0 , w = 0 . Equation then becomes 



dh dh 



- + - 



dx^ dy^ 



= u 



^ pc'\ dt 



V «• y 



dx 



(3.9) 



a = kj pc, the thermal diffusivity of the fluid. So the governing equation for the parallel 
plates flow is 



dh dh u dt 



dyr dy^ a dx 



(3.10) 



Now, to non-dimensionalize the problem, assume that 



R-dlZL V+-Z 
h-h 






/n-r. h' V' RePr ' 2 



The governing equation becomes 



1 d^O ^ d^e u* de 



Pe 



■ + 4 



Dh dx" 



dy^ 



2 



(3.11) 



13 



4. Flow Between the Parallel Plates, Constant Heat Flux Case 



The governing equation for the parallel plates flow is (3.10) 



d^t d^t u dt 



■ + ■ 



dx'^ dy a dx 



Now, to non-dimensionalize the problem, assume that 



+ T + “ + 2x1 + 3^, +2 

q^Djk h r’ RePr ’ 2^ 






The governing equation becomes 



, d^e u* de 
- + 4- 



^ 4 , 5 .*^ 2 &• 



(3.12) 



5. Effects of Wall Conduction, Flow Inside the Circular Tube, 



Constant Wall Temperature 



The governing equation for the circular duct flow is (3.4) 



d^t 1 dt 8^t u dt 



- + + 



dr^ r dr dx^ a dx 



Now, to non-dimensionalize the problem, assume that 

tQ-t ^ r ^ u 4 . x/r, 



e=- 



h-te 



+ r + w 

r = — , u = — , X = 

rn V RePr 



+ - •'"O , u* =2(l-r* ) 



The governing equation becomes 



l-r+ 






d^0 1 d& 1 d^0 



■ + ■ 



dx^ dr*' dr* Pe^ dx* 



(3.13) 
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B. BOUNDARY CONDITIONS 



1. Flow Inside the Circular Tube, Constant Wall Temperature 
Case 




Figure 3. The boundary conditions of circular tube, constant wall temperature flow 



The governing equation is (3.5) 

Aee d^e i ee i d'^e 

^ ^ dr* r* dr* Pe^ dx* 

where the parameters of the equation are defined above. The boundary conditions for this 
case is: 



at r^=0, 


-^ = 0 
dr* 


(Symmetry boundary condition) 


(3.14a) 


at r'^=l, 


0 = 0 


for 

for 


+ + 
A IV 
O O 


(3.14b) 


as X* - 00 , 


t ->tey 


(9->l 




(3.14c) 


as X* ->+ 00 , 




6>->0 




(3.14d) 
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2. Flow Inside the Circular Tube, Constant Heat Flux Case 



Insulated wall r 

////// \ 


U 


_J 


u 


1 1 




r* A q = Qo = constant 

X=-00 ^ no x = +oo 


^ = 0 ^ ^ e-^-2x* 

Constant beat flux 


////// 
nsulated wall 


^ i 




, 


> 





Figure 4. The boundary conditions of circular tube, constant heat flux flow 



The governing equation is (3.6) 

L +2)36^ d'^0 1 36 1 3'^ 6 

^ ^ dr* r* dr* Pe^ 3x* 

where the parameters of the equation are defined above. The boundary conditions for this 
case is: 

at r"^=0, 

at r^=l, 

as X* — > — CO , 
as ^ + 00 , 



ae 

dr* 



= 0 



ae 1 

8r* 2 

ae 



.ar* 



•=o 



(Symmetry boundary condition) 



for X* >0 

X* <0 



t , 6 — > 0 



6 — ^ — '2.x + 



+4 ,2 

7 r* r* 

) 

/IQ 5 9 



(3.15a) 

(3.15b) 

(3.15c) 

(3.15d) 
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3. Flow Between the Parallel Plates, Constant Wall Temperature 



Case 



y Insulated wall 

yf yf y y y Tv«a = To = constant 




, — ► y*t , 


X->+CO 


^ 1 y u(y^) 

/ _ 


^ ^->0 


^ ^ Constant wall temperature 

Insulated wall 





Figure 5. The boundary conditions of parallel plates, constant wall temperature flow 



The governing equation is (3.1 1) 
1 






+ 4- 









u* de 

2 dx^ 



where the parameters of the equation are defined above. The boundary conditions for this 
case is: 



o' 

II 


^ = 0 


(Symmetry boundary condition) 


(3.16a) 




f(9 = 0 


for 

for 






aty^=l, 


II 

o 


X* >Q 

X* <0 


(3.16b) 


as -oo 5 




e-^\ 




(3.16c) 


as + 00 , 


^ ^05 


9-*Q 




(3.16d) 
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4. Flow Between the Parallel Plates, Constant Heat Flux Case 



y Insubted wall r 

////// 1 


U 


. i 1 


u 




y* A q = constant 

x-F-oo ^ T ^ ^ x = +oo 


y u(y") x’* ^ Q^^2x* 

y Constant heat flux 


////// 
Insulated wall 




0 ft 







Figure 6. The boundary conditions of parallel plates, constant heat flux flow 



The governing equation is (3.12) 

L +2\d0 d^e 1 dO 1 d^e 

^ ^ dr^ dr^ Pe^ dx^ 

where the parameters of the equation are defined above. The boundary conditions for this 
case is: 



+ 30 

at y^=0, — = 0 (Symmetry boundary condition) 

dr* 



aty^=l. 



de _ 1 

” 4 

de 



dr* 



• = 0 



for X* >0 
x"<0 



2LS — > —00 j t — > ^ 0 — > 0 



HS — > +00 , t — > /q j 0 — — + 



39 3 

- + ■ 






1120 16 



-y 



(3.17a) 



(3.17b) 



(3.17c) 

(3.17d) 
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5. Effects of Wall Conduction, Flow Inside the Circular Tube, 



Constant Wall Temperature 




Figure 7. The boundary conditions of circular tube, constant wall temperature including the circular wall 
The governing equation is (3.5) 






d^e 1 dd 1 d^e 



dx^ ^ dr^ ' Pe^ dx^" 



The boundary conditions for this case is: 

dO 



at r"^=0, 



dr* 



= 0 (Symmetry boundary condition) 



(3.18a) 



at r^=l, 



at r^=l+<5r , 



^watl ^ fluid 

^^wall _ fluid 



dr* dr* 
' d9 



dr* 

9 = 0 



= 0 



for X* <0 
x*>0 



(3.18b) 



(3.18c) 



where c5r is — , r* is the thickness of the wall, and ro is the inner radius of the tube. 



as 


x"" -> -00 , 




e 


->i 


as 


->+00, 




6 


->o 



(3.18d) 

(3.18e) 
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IV. NUMERICAL REPRESENTATION 



A. FLOW INSIDE THE CIRCULAR TUBE, CONSTANT WALL 
TEMPERATURE CASE 



1. The Governing Equation 



The governing equation is (3.5) 



where 




dd 

dx* 




1 80 1 d^0 





^ x/rp 

RePr 



The finite difference approximations for the first and second derivatives with error of 
order h^ are (figure 2.) 



dx\,n 2K 



(4.1a) 



de 

Sr* (ij) 



2h 



+ o{hl) 



(4.1b) 



d-9 

Sx*~ (,j) 



hi 



■oihl) 



(4.1c) 



21 



(4.1d) 



d^e 

(U) 






When these approximations are plugged into the governing equation (3.5), the equation 
becomes 



(4.2) 



1 




+ 




1 


f ff ^ 

^(^7+1) “ J-1) 


-(’->'(1.)] 




Pel 


1 ti J 


J 




2A, J 


1 2 ^ J 



After simplifying the expression (4.2), becomes 









Pelhl 



yl.j) 

2h. 



J_ 

2K 



+ 0 , 






1 



^ + J_ 

Pelh^ 2A. 2A. 



+ 6( 






yK 












OJ~l) 



hi' 



'^Ky<uj: 



2 2 
PeW.^'hl 



(4.3) 



2. The Boundary Conditions 

The boundary conditions are: 



a) At the Center 

+ 30 

At r~^=0. = 0 (Symmetry boundary condition) 



1 36 

When this boundary condition is applied to the governing equation, the element 

r* dr* 

becomes ^ . After applying the L'Hospital's rule, the governing equation for the given 
boundary condition becomes 
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(4.4) 



^-2-^ ^ 

fix* ^ ax:*' 

Since the temperature distribution at the center on the radial direction is symmetric, 



After applying the finite difference method, the 9^^ at the boundary becomes 



9, 



m.j) 






Pp^ / j 2 2h 






+ ■ 



Pp^ 2 h 






2 ,,2 + ^2 

y 



p^dK 






(4.5) 



b) At the Wall 



atr^=L 



^(U)=0 

39 



lar* 



= 0 



for X* > 0 

fo’' x*<0 



(4.6) 



FOR x^<0 



x^<0, r^=l, — = 0 

ar* 

The governing equation for the given boundary condition becomes 

3^9 1 3^9 ^ 

— r + — ; r = 0 

ar* Pel 5x* 



Since the wall is insulated. 



39 



3r 



+ - 0 > ^(1 j+i) “ ^{ij- 



at the boundary. 



(4.7) 



After applying the finite difference approximations, the at the boundary becomes 
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^0+1.;) 






J 






Pelhlj 



+ 9. 



^ 2 ^ 



(>.;-i) 



v^;y 



2 ^ 
Pelh; ^ h; 



(4.8) 



c) As X -^ - 00 

asx"-^-», (4.9) 

d) As X “^ + 00 



3.S — > +00 5 0^. ^ 0 



(4.10) 



3. The Bulk Mean Temperature 

The mean temperature at a position “i” is (2.8) 

2 

=-^ K,r)t(x,rfdr 
^^0 0 

'b 

where V is the average velocity and to is the radius of the tube. When \t^ul7trdr is 

0 

subtracted from both sides and dimensionless parameter r^ is substituted inside the 
equation, the equation becomes 

{tm-h)= \{t-t^)u*2r^dr^ 

0 

Let’s define the bulk mean temperature, 9^ = ^ , (4.11 a) 



24 



where te is the inlet temperature (as x — > -oo ) and to is the constant wall 
temperature.(Fluid exit temperature as x ^ +oo ) 

Then the bulk mean temperature 9^ becomes 

9 ^ =4jfl-r^^V + (4.11b) 

>n(x+) JV 

9 is found out from the numerical calculation. After that is calculated by solving the 
equation (4.1 1) numerically. 



4. The Nusselt Number 

The Nusselt number is (2.3) 




dO 

When the heat transfer equations = h{Tg -T^) and = -k — 

dr 

dimensionless parameter r^ are plugged in (2.3), the equation becomes 



, and the 



r=R 



Nu = - 



2 69 



9„ dr'" 



r+=l 



(4.12) 



The finite difference approximation for the equation (4.12) using the backward 
difference expression with error of order h becomes 



= — 
(0 n 



m(i) 



2h.. 



(4.13) 



where “j” is the value which corresponds to r^=l . 
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5. Heat Absorption Upto the Entrance 

Assume a differential volume, dV = iTirdr between x = -oo and x = 0 . The 
energy stored in this control volume is 

dE = CpAtdm (4.14) 

where At is the temperature difference between entrance temperature and the steady state 
temperature after the heating starts, At = - 1^ and dm is the mass of the control 

volume, dm = pdV .The energy absorbed in the control volume becomes 

0 ’’o 

^absorbed = p J |(^ “ te)rdrdx . (4.15) 

-CO 0 



0 ''o 

When 2npCp | \tjdrdx is subtracted fi'om both sides and (te-to) is divided to both 

—00 0 

sides, the equation (4.15) becomes 

p 0 '’o 

^absorbed ^ f \{6 -\)rdrdx (4.16) 

ite-to)pCp -io 



t — t xiv T 

where 9 = As the non-dimensional parameters of and r* = — and 

RePr r„ 



te-to 



their first derivatives are substituted into the equation (4.16), the equation becomes 



p 0 1 

. . = _2nPe f \{9 -ly^dr-^dx-^ 



(4.17) 
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B. FLOW INSIDE THE CIRCULAR TUBE, CONSTANT HEAT 
FLUX CASE 



1. The Governing Equation 



The governing equation is (3.6) 



where 




. 2 \ 



89 



dx* 




1 80 1 8^0 
^r^8r^^ Pel 






ql d/k 




+ _ 2x1 d 
Re^ Pr ■ 



The finite difference approximations for the first and second derivatives with the error of 

'y 

order h are (4. 1 a,b,c,d) (figure 2.) 



80 

^ (ij) 



2h 



+ o{hl) 



80 

dr* (ij) 



2h.. 



+ o{hl) 



8^0 

dx (iJ) 






h 



2 

X 



+ o(h^J 



8^0 

dr*' o.» 






+ o{hl) 
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When these approximations are plugged into the governing equation (3.6), the equation 
becomes 



^(»+U) ^(/-U) ^ ^^(Ij) 1 


V 




^ 1 








} 


J 




1 J 


1 2/,, J 



= 0 



(4.18) 



After simplifying the expression (4.18), becomes 






o I 1 , yl.n 1 

<'*'•^>1 Pe ^ t ] 2K 2 / 1 , 






(i-W) 



Pey,] 



Ih 2h 



2V«.. 



+e, 






>J)J 



hi 



IK 






2 2 
Pelh]^ hi 



(4.19) 



2. The Boundary Conditions 

The boundary conditions are: 

a) At the Center 

+ 30 

At r~^=0. = 0 (Symmetry boundary condition) 

dr* 

1 dG 

When this boundary condition is applied to the governing equation, the element 

r* dr 

becomes ^ . After applying the L'Hospital's rule, the governing equation for the given 
boundary condition becomes 



^ = 2— + — — 
dx* ~ dr*‘ Pe^ dx*' 



(4.20) 
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Since the temperature distribution at the center on the radial direction is symmetric, 






After applying the finite difference method, the 6^^^ at the boundary becomes 



9, 









1 1 



Pe'h^ 2h 






1 1 



Pp^ 2h 

\repfix J 



+^, 









2 2 
+ ■ 



Pe\hl h; 



(4.21) 



b) At the Wall 



at r^=L 



FOR >0, 



d9 

dr* 



2 



k^ = q'o, 9{x*,r*) = 
or 



for 


>0 


for 


+ 

A 

O 


K ~ ^(x,r) 




q'odlk 


■M 

r 



r 

r 



(4.22) 



u • u dO \ 

so the equation becomes = — . 

dr* 2 



The finite difference approximation for this equation is 






2K. 



When the boundary condition is applied to the governing equation (3.6), the equation 
becomes 

d^9 1 1 d-9 



dr*" 2r* Pe- dx*" 



= 0 



(4.23) 
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As the finite difference approximations are substituted in the equation above, the equation 
becomes 









1 ^ 

yPelhlj 






1 ^ 






+ 9 > 



p^dK K 



FOR X <0 






1 1 



V‘y J ‘y 



h„ 2 



x^<0, r^=l, 



39 

dr* 



= 0 



The governing equation for the given boundary condition becomes 



3^0 1 3^9 ^ 

+ : r- = 0 






dO 

Since the wall is insulated, = 0 , = ^(,,y_i) at the boundary. 



After applying the finite difference approximations to (4.25), the 9,^ 






becomes 



9 . 



0+1.;) 



P(iJ) ~ ■ 



1 ^ 



+ ^o-l. 



J) 



1 ^ 

s^Pelhlj 



+ 9 , 



^2^ 






\J^y J 



2 2 
+ 



Pelh] h; 



c) As X - 00 



as X — ^ — CO , 9^^ jy — ^ 0 



(4.24) 



(4.25) 



at the boundary 



(4.26) 



(4.27) 



30 



d) As X -> + 00 



3S X — ^ +CX) j 



^(ij) ■*■ 



7 



• + ■ 



48 8 







(4.28) 



3 . The Bulk Mean Temperature 

The mean temperature at a position “i” is (2.8) 



2 

^ jutrdr 
^'0 0 



fo 

where V is the average velocity and ro is the radius of the tube. When \tgUl 7 rrdr is 

0 

subtracted from both sides and dimensionless parameter r^ is substituted inside the 
equation, the equation becomes 

(4.29) 

0 

Let’s define the bulk mean temperatoe, 



" J/7 ’ 



(4.30a) 



where te is the inlet temperature (as x -> -co ). 
Then the bulk mean temperature 9^ becomes 

0 



(4.30b) 
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0 is found out from the numerical calculation. After that 6^ is calculated by solving the 
equation above numerically. 



4. The Nusselt Number 



The Nusselt number: iVw = — . 

k 

The heat transfer from the wall: q\ = h{t . » ~ • 

The temperature at the position of (x’^,r'^=l): t . , =^e 

The average temperature at (x^): =t^ 

When these equations are substituted into the Nusselt number equation, Nusselt 
number becomes: 



Nu = 



1 

6 -6 
^ m ^ 



(4.31) 



5. Heat Absorption Upto the Entrance 

Assume a differential volume, dV = Ijo-dr between x = -oo and x = 0 . The 
energy stored in this control volume is (4.14) 

dE = Cpbddm 

where At is the temperature difference between entrance temperature and the steady state 
temperature after the heating starts. At = and dm is the mass of the control 

volume, dm = pdV .The energy absorbed in the control volume becomes 
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0 ''o 

^absorbed = p j \{t ~ te)rdrdx . 
—00 0 



(4.32) 



When both sides are divided by qgdjk , the equation (4.32) becomes 



p 0 ''o 

. = -Ak j \6rdrdx 



-00 0 



(4.33) 



t — t X r T 

where 6 = — ^ . As the non-dimensional parameters of =_L_°_ and y = — and 

qyjk RePr 

their first derivatives are substituted into the equation (4.33), the equation becomes 



^absorbed ^ | + + 



^PCpKlk 



(4.34) 



-00 0 



C. FLOW BETWEEN THE PARALLEL PLATES, CONSTANT 



WALL TEMPERATURE CASE 



1. The Governing Equation 

The governing equation is defined in equation (3.1 1) as 



Pe 



1 d^e , d^e y de 

— - + 4- 



Dh dx~^ dy 



y 2 dx^ 



where 



^ t(.-t + y + 2x1 Du ^ 

0 = - , y =T> ^ D^=Ah. 



to-te 



h 



RePr 



The finite difference approximations for the first and second derivatives with error of 



order m are (4.1a,b,c,d) (figure 2.) 
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ax^ IK 



dO 

dr* (ij) 



2h„ 



+ o(K,) 



d~e 

2 

dx (,j) 



^('+1.7) ^0-hj) 



- 20 . 



(iJ) 



h 



2 

X 



+ o{h;) 



d^O _ ^(i,;+l) + ^(.,7-1) 2^(/,7) . /7 2 x 

T Ti ^ 0 (hy) 

dr* oj) ^7- 

When these approximations and u"^ are plugged into the governing equation (3.11), 0^^ 
becomes 






^(■•7) “ ■ 






Pel hi Sh, 



Sh, 



+ 0 . 



(>-1.7) 



Pel.hl 



, 3 
8/!. %h. 



\ 


f4^ 






f4l 


/ 






>-l) 





(4.35) 



2 

Pel hi ^ hi 



2 . The Boundary Conditions 

The boundary conditions are: 

a) At the Center 

dd 

At y"^=0, j=l, = 0 (Symmetry boundary condition) 

dy* 
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When this boundary condition is applied to the governing equation, the governing 
equation becomes 



1 d'^e , d^e 3 dG 
- + 4- 



(4.36) 



Pe\ cbc*' dy*" 4 dx* 

Since the temperature distribution at the center on the radial direction is symmetric, 

at the center. 

After applying the finite difference method to (4.36), the at the boundary becomes 



G, 



(/+1.1) 



^(M) “ ' 



1 



^Pelhl 






1 3 

+ - 



^Pelh; ZK 







+ ^0'.2) 




j 


\^y 7 



2 ^ 

Pelhl ^ hi 



(4.37) 



aty"=l. 



b) At the Wall 

for X* >0 

dy* for X* <0 



dG 



= 0 



(4.38) 



FOR x!<0 



x"<0, y"=l. 



dG 

dy^ 



= 0 



The governing equation for the given boundary condition becomes 

1 d-G , d-G ^ 

• + 4 r = 0 



Pei dx*" dy*" 



(4.39) 
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dO 

Since the wall is insulated, — ^ = 0 , .j, at the boundary. 

dy* 

After applying the finite difference approximations, the .. at the boundary becomes 






{MJ) 



^OJ) ~ ' 



1 ^ 



+0, 



0-i.» 



1 ^ 



+ 0 , 



^8^ 



2 8 
+ - 



Pelh^, h; 









(4.40) 



e) As X ^ - 00 



as X — ^ — CO , 0 ^ 1 



(4.41) 



f) As X "^ + 00 



as X* ->+oo, 0(ij) 



(4.42) 



3. The Bulk Mean Temperature 

The mean temperature is defined in equation (2.9) as 

c=^W. 

When the dimensionless parameters u^, y^ are substituted in (2.9), the tm becomes 



= \u^tdy^ 

0 

where V is the average velocity and h is the distance fi'om the centerline to the plate. 
1 

When = \t^u^dy^ is subtracted from both sides, the equation becomes 

O' 
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0 

Let’s define the bulk mean temperature, 

(4.43a) 

where te is the inlet temperature (as x -oo ) and to is the constant wall 
temperature.(Fluid exit temperature as x ^ +co ) 

Then the bulk mean temperature 6^ becomes 

^ 0 

9 is taken fi'om the numerical calculation. After that 0^ is calculated by solving the 
equation (4.43) numerically. 



4. The Nusselt Number 

The Nusselt number is defined as (2.3) 



Nu = 



k 



dO 

When the heat transfer equations q = h(T -T^), ql = k — 

dy 



, and the dimensionless 



y=h 



parameter r^ are plugged in (2.3), the equation becomes 



Nu = - 



4 dO 

e.dyy 



(4.44) 



/=! 
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The finite difference approximation for the equation (4.44) using the backward difference 
expression with error of order h becomes 



Nu^ = 






mi 



2h„ 



(4.45) 



5. Heat Absorption Upto the Entrance 

Assume a differential volume, dV = bdxdy between x = -oo and x = 0 . The 
energy stored in this control volume is (4.14) 

dE = CpEtdm 

where Ar is the temperature difference between entrance temperature and the steady state 
temperature after the heating starts, At = y-^ - and dm is the mass of the control 

volume, dm - pdV .The energy absorbed in the control volume becomes 

^absorbed ~ J (4-46) 

-00 -/2 

0 h 

When 2bpCp J \t^dydx is subtracted fi’om both sides and is divided to both 

-00 0 

sides, the equation (4.46) becomes 

_ Eabsorbed ^ 2 ] \d ~ \)dydx (4.47) 

bpCp{t,-Q 
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where 0 = . 

iK-to) 



As the non-dimensional parameters of 



2x!D, 

RePr 



and 



T 

= — and their first derivatives are substituted into the equation (4.47), the equation 
becomes 



F 0 1 

= -APe j \{G-\)dy*dx^ 



h bpcpit^-t^) 



(4.48) 



-00 0 



D. FLOW BETWEEN THE PARALLEL PLATES, CONSTANT 
HEAT FLUX CASE 



1. The Governing Equation 



The governing equation is (3.12) 



1 d^e , d^e u* 30 
+ 4- 



dx* dy^ 2 



where 



t.-t 
0=^ 



^_y 



llDJk’ ^ h’ ^ RePr 



The finite difference approximations for the first and second derivatives are (4.1a,b,c,d) 



30 

3x* (ij) 



2h, 



30 

dr* (ij) 



2h ^ 
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d^e 



(ij) 



h'x 



d^e 
dr^\i.n 









+ o(/z‘) 



When these approximations and are plugged into the governing equation, the bulk 
temperature at the position (i J) is wnritten as 






^(ij) ~ ■ 






yPelhl Sh, 









_! ^yi.j ) , 3 

^Pelh] 8h^ 8K 



+e, 



(i,y+i) 









+ 6 . 






^_4_^ 

yK, 



(4.49) 



2 8 
PeW^hl 



2. The Boundary Conditions 

The boundary conditions are: 

a) At the Center 

do 

At y^=0, j=l, = 0 (Symmetry boundary condition) 

dy* 

When this boundary condition is applied to the governing equation (3.12), the governing 
equation becomes 

1 d-e . d‘6 3 d9 .. ... 

— - + 4 — ^ = (4.50) 

Pe\ 4 

Since the temperature distribution at the center on the radial direction is symmetric, 
*e center. 
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After applying the finite difference method, the at the boundary becomes 






1 3 



^(M) = 



^Pelhl 8 / 2 ,^ 






1 ) 



Pelhl Bh 



+ - 






y 






8 



^4^; K 






(4.51) 



b) At the Wall 



aty =1, = < 



0 

_J_ 

4 



for x"" < 0 



for X* >0 



FOR x^<0 



x^<0, y"=l, 



de 

dy* 



= 0 



The governing equation for the given boundary condition becomes 

1 d‘0 . d-e - 

■ + 4 r = 0 



dx* dy* 



(4.52) 



dd 

Since the wall is insulated, = 0 , 0.^ = 0.^ _,s at the boundary. 

dy* 

After applying the finite difference approximations, the 0,^ ^ at the boundary becomes 






^(>j) ~ 



1 ^ 
^Pelh:^ 



+ 0 , 



(i-lj) 



1 ^ 

y Pelhl J 



+ 0 , 












2 8 



(4.53) 



Pelhl hi 
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FOR x" > 0 



x*>0, y‘'=l, 



60 ^ 1 
" 4 



The governing equation for the given boundary condition becomes 

1 6-0 , 6^0 ^ 

■ + 4 r = 0 



6x*‘ 6y* 



(4.54) 



6y* 4 ’ 2h 4 



= ^( 0 - 1 ) ■ Y boundary. 



After applying the finite difference approximations, the at the boundary becomes 






^0 j) ~ ■ 



1 ^ 






1 ^ 
KP^lhlj 









+ • 






(4.55) 



c) As X - 00 



as -^-oo , 0^^ — > 0 



(4.56) 



d) As X -> + 00 



as -> + 00 , 0(^. = -2x^ + 



39 3 

• + 



1120 16 



.4 

V 



(4.57) 



3. The Bulk Mean Temperature 



The mean temperature is (2.9) 

1 ^ 

tm =-^Wdy. 



Vh, 
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When the dimensionless parameters u"^, are substituted in (2.9), the tm becomes 



t. = 

0 

where V is the average velocity and h is the distance from the centerline to the plate. 

1 ^ 

When is subtracted from both sides, the equation becomes 

0 

0 

Let’s define the bulk mean temperature, 



0 ^ = 






(4.58a) 



where te is the inlet temperature (as x — > -oo ). Then the bulk mean temperature 



becomes 



^ 0 

6 is taken from the numerical calculation. After that 6^ 
equation (4.58) numerically. 



(4.58b) 

is calculated by solving the 



4. The Nusselt Number 

hD 

The Nusselt number; Nu = — - . 

k 

The heat transfer from the wall : ql = h{t , , ~K) ■ 

(JT ,r —1) 

The temperature at the position of (x^,r^=l): D^Jk 
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The average temperature at (x^: 

When these equations are substituted into the Nusselt number equation, Nusselt number 
is calculated as 



A^w = 



1 

e -e 

^ m 



(4.59) 



5. Heat Absorption Upto the Entrance 

Assume a differential volume, dV = bdxdy between x = -oo and x = 0 . The 
energy stored in this control volume is (4.14) 

dE = CpAtdm 

where At is the temperature difference between entrance temperature and the steady state 
temperature after the heating starts, At = and dm is the mass of the control 

volume, dm = pdV .The energy absorbed in the control volume becomes 



0 +/? 

^absorbed = J " t,)hdydx . 

-OO -h 



(4.60) 



II / ^ 

When both sides are divided \sy q^D^jk , the equation (4.60) becomes 



r 0 /j 

-Elabsprbed_^_^ I 

-coO 



(4.61) 



where 0 = ^ 



2^ / 7~) f* 

As the non-dimensional parameters of x"^ = — - — - and = — and 
qoDjk RePr 



their first derivatives are substituted into the equation (4.61), the equation becomes 



/r 0 1 

^ . f b ^ped . ^ _i^pg J \Sdy*dx^ 

iJjh ja _„o 



(4.62) 



44 



E. THE EFFECTS OF WALL CONDUCTION, FLOW INSIDE THE 



CIRCULAR TUBE, CONSTANT WALL TEMPERATURE CASE 



1. The Governing Equation 



The governing equation is (3.5) 



where 



r, ee 

1 -r 

y 



d^0 



dx^ 



■ + ■ 



1 de 1 






d^e 



dr^ 



dr^ Pe 



D dx^ 



e = 



h-K ’ 




RePr 



The finite difference approximations for the first and second derivatives with error of 
order h are (4.1a,b,c,d) (figure 2.) 



&*,u) -i-K 



de 

dr\u) 



2h„ 



+ o{hj.) 



d^e 



9. 



O+l.j) 






-2^0.) 



h 



2 

X 



■^o{hl) 
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d-e 

dr*' ii.j) 






^o{hl) 



When these approximations art plugged into the governing equation (3.5), the equation 
becomes 









-26. 



VJ) 






>'0.y) 



^(1 J-n) ~ ^0J-\) 



2K 






ylj) 









2K 



= 0 



(4.63) 



After simplifying the expression (4.63), the equation becomes 









Pelh] 



yjj) 

2h 



2h 



+ 0. 






Pelh] 



ylj) 



1 



2h 2h 



+ ^r 



(t.y+i) 



1 



2V(../ 



+e, 



(i.y-i) 



h] 



^Ky«.j 



(4.64) 



2 2 
PeW. ^ 



2. The Boundary Conditions 

The boundary conditions are: 



a) At the Center 

4 - 30 

At r"^=0, — ^ = 0 (Symmetry boundary condition) 



When this boundary condition is applied to the governing equation, the element — - — - 

r dr 

becomes ^ . After applying the L'Hospital's rule, the governing equation for the given 
boundary condition becomes 



dO d‘0 1 d-e 

2* ~h ' 



dr*" Pe' dx*" 



(4.65) 



Since the temperature distribution at the center on the radial direction is symmetric. 
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After applying the finite difference method, the s at the boundary becomes 



e, 



('+i.y) 



- ' 



Pe>^ 2h , 






J) 



1 1 

• + ■ 



PeW 2K 



Pelh] ■" h] 





^ 4 ^ 


+%7+l) 




y 


\^y ) 



at r^=l, 



b) At the Wall 



^wall ^ fluid 

^^wall _ fluid 



[ dr^ dr^ 

The heat fluxes through the fluid and the wall are equal to each other at r^ 



i, ^^solid _ K fluid 

solid : “ 



dr^ dr* 

Using the 2"‘* order forward and backward difference method, the 
becomes 





- k 


^^OJ) ~ ^(> J+2) ^ 






^»y. J 



If hys=hyf_ the bulk temperature at the fluid- wall boundary becomes 









%/)- 

1 + -^ 

V ^sj 



w kf 

Oti + ~ — dit 



'i‘J+2) 






J 



(4.66) 



= 1 . 



equality 



(4.67) 



(4.68) 
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c) On the Outer Edge of the Wall 



FOR x^<0 



For x"^<0 and r'^ = 1 + , the outer surface of the solid wall is insulated where 



dy = 



solid _ wall _ thickness 
tube inner radius 



The governing equation for the given boundary condition becomes 



d'^e 1 d^e ^ 

+ r r = 0 



dr*' Pe^ dx*' 



(4.69) 



dO 

Since the wall is insulated, = 0 , at the boundary. 



After applying the finite difference approximations, the at the boundary becomes 






0*U) 



~ ■ 



1 ^ 






1 ^ 






^2^ 






2 2 
+ 



Pelh^. h; 



(4.70) 



FORx^>l 



^ihj) - ^ 



(4.71) 



d) As X ^ - 00 

as X* -> - 00 , G^ij^ 1 (4.72) 

e) As X ^ + 00 
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3^S X — ^ -hoo 5 



(4.73) 



0. 



(<j) 



^0 



3. The Bulk Mean Temperature 

The mean temperature at a position “i” is (2.8) 

2 

0 



ro 

where V is the average velocity and ro is the radius of the tube. When \t^u2mdr is 

0 

subtracted from both sides and dimensionless parameter r"^ is substituted inside the 
equation, the equation becomes 

0 

Let’s define the bulk mean temperature, ^ ^ where te is the inlet temperature 

V^c ~^o) 

(as X -> -00 ) and to is the constant wall temperatm’e.(Fluid exit temperature as x — > +oo ) 
Then the bulk mean temperature becomes 



6 . 



m{i) 




(4.74) 



0 is found out from the numerical calculation. After that 6^^ is calculated by solving the 
equation (4.74) numerically. 
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4. The Nusselt Number 



The Nusselt number is defined as (2.3) 



Nu= — . 
k 



II II 

When the heat transfer equations = h(t^ -t^) = -k 



dimensionless parameter r^ are plugged in (2.3), the equation becomes 



dr 



, and the 



r=R 



Nu=- 



2 80 



0„ dr* 



r+=l 



(4.75) 



The finite difference approximation for the equation (4.75) using the backward 
difference expression with error of order h^ becomes 



^m{i) 



2h 



(4.76) 



where “j” is the value which corresponds to r^=l. 



5. Heat Absorption Upto the Entrance 

Assume a differential volume, dV = 2m-dr between x = -oo and x = 0 . The 
energy stored in this control volume is (4.14) 

dE = Cplstdm 

where At is the temperature difference between entrance temperature and the steady state 
temperature after the heating starts. A/ = and dm is the mass of the control 

volume, dm = pdV .The energy absorbed in the control volume becomes 



50 



0 '’o 

E absorbed = p J |(^ " te)rdrdx . 

—00 0 



(4.77) 



0 '■(> 

When iTtpCp j \t^rdrdx is subtracted from both sides and (te-to) is divided to both 

— 00 0 



sides, the equation (4.77) becomes 



f 0 ''o 

— sh^orb^ ^ 2^ I J(^ _ \ydrdx 

^e~^t^P^p -« 0 



(4.78) 



t ^ t x! r Y 

where 6 As the non-dimensional parameters of x'*' = ■■■ ° and = — and 

RePr r. 



te-to 



their first derivatives are substituted into the equation (4.16), the equation becomes 



p 0 1 

^absorbed ^ f UQ -\Y dv^ dx^ 



(4.79) 
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V. NUMERICAL RESULTS AND DISCUSSION 



A. BULK MEAN TEMPERATURE VS. X 

The bulk mean temperature graphs of the fully developed laminar flow inside the 
circular tube and parallel plates for constant wall temperature and heat flux cases were 
plotted in figures 8-11. The Peclet numbers 0.5, 1, 2, 3, 5, 10, 50, and 100 were used for 
the plots to see how the Peclet number, Pe=RePr, affects the axial conduction in the flow. 
As it is seen from the figures, when the Peclet number gets smaller, the axial conduction 
effects of the fluid becomes very important and the temperature of the fluid upto the 
entrance increases. The axial conduction of the fluid has to be considered for these cases 
and must not be neglected for the Peclet numbers less than about 10. When the Peclet 
number gets larger, the effects of the axial conduction decrease. As it is seen from the 
figures 8-11, for the Peclet numbers 50 and 100, the bulk mean temperature of the fluid at 
the position x‘"=0 is almost equal to the entrance temperature, which means that there is 
almost no heat absorption or temperature increase upto the entrance. The axial conduction 
of the fluids can be neglected for these high Peclet numbers. 

1. Flow Inside the Circular Tube, Constant Wall Temperature 

Case 

The graphs of the bulk mean temperature vs. x for the circular tube, constant wall 
temperature flow was plotted for various Peclet numbers in figure 8. a and 8.b where theta 
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bulk mean is 






{tm -to) 

iU-t.) 



and the nondimensional axial distance x* is x* 



2x1 d 
RePr 



The validity of the figure 8.b has been checked by comparison with Ref.5. 



CIRCULAR TUBE CONSTANT WALL TEMP. FLOW, THETA BULK MEAN 




Figure 8.a The axial bulk mean temperature distribution for circular tube constant wall temperature flow 



CIRCULAR TUBE CONSTANTWALL TEMP. FLOW. THETA BULK MEAN 




Figure 8.b The enlarged view of figure 8.a 
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2. Flow Inside the Circular Tube, Constant Heat Flux Case 



The graphs of the bulk mean temperature vs. x for the circular tube, constant heat 
flux flow was plotted for various Peclet numbers in figure 9.a and 9.b where theta bulk 

mean is 9„ = and the nondimensional axial distance x^ is x* = ^ _ 

q^dlk Re^Pr 



CIRCULAR TUBE CONSTANT HEAT FLUX FLOW. THETA BULK MEAN 




Figure 9.a The axial bulk mean temperature distribution for circular tube constant heat flux flow 

CIRCULAR TUBE CONSTANT HEAT FLUX FLOW. THETA BULK MEAN 




Figure 9.b The enlarged view of figure 9.a 
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3. Flow Between the Parallel Plates, Constant Wall Temperature 
Case 

The graphs of the bulk mean temperature vs. x"" for the parallel plates, constant wall 
temperature flow was plotted for various Peclet numbers in figure lO.a and lO.b where 

theta bulk mean is ^ ~ ) and the nondimensional axial distance x"^ is x* = • 

” RePr 



PARALLEL PLATES CONSTANT WALL TEMP. FLOW, THETA BULK MEAN 




Figure lO.a The axial bulk mean temperature distribution for parallel plates constant wall temperature flow 



PARALLEL PLATES CONSTANT WALL TEMP. FLOW. THETA BULK MEAN 




Figure lO.b The enlarged view of figure lO.a 
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4. Flow Between the Parallel Plates, Constant Heat Flux Case 



The graphs of the bulk mean temperature vs. x for the parallel plates, constant heat 
flux flow was plotted for various Peclet numbers in figure 1 1 .a and 1 1 .b where theta bulk 

mean is ^ and the nondimensional axial distance x"^ is 

qoDjk RePr 

PARALLEL PLATES CONSTANT HEAT FLUX FLOW. THETA BULK MEAN 




Figure ll.a The axial bulk mean temperature distribution for parallel plates constant heat flux flow 

PARALLEL PLATES CONSTANT HEAT FLUX FLOW, THETA BULK MEAN 




Figure ll.b The enlarged view of figure 1 l.a 
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B. NUSSELT NUMBER VS. X 



The nusselt number graph of the fully developed laminar flow inside the circular 
tube and parallel plates for constant wall temperature and heat flux cases were plotted in 
figures 12-15. The Peclet numbers 0.5, 1, 2, 3, 5, 10, 50, and 100 were used for the plots 
to see how the Peclet number, Pe=RePr, affects the nusselt number distribution in the 
flow. 

1. Flow Inside the Circular Tube, Constant Wall Temperature 

Case 

The graph of the nusselt number vs. x for the circular tube, constant wall 
temperature flow was plotted for various Peclet numbers in figure 12 where x"" is 
+ 2x1 d 

X = . 

RePr 

CIRCULAR TUBE CONSTANT WALL TEMP. FLOW, NUSSELT NUMBER 
30 

25 

20 
V 

n 
E 

2 15 



5 

10 '^ 10 ’^ 10 ° io’ 

X 

Figure 12. The nusselt number distribution for circular tube constant wall temperature flow 



CIRCULAR TUBE CONSTANT WALL TEMP. FLOW, NUSSELT NUMBER 
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The nusselt numbers approach to the values in those of table 1 when x gets larger as seen 



in figure 12. 



Pe 


0.5 1 


2 


3 


5 


10 


50 100 


Nu. no 


4.0971 4.0276 


3.9224 


3.8506 


3.767 


3.6948 


3.6586 3.6572 


Table 1. T 


le nusselt number for various Peclet numbers for circular tube constant wall temperature case as 



X goes to infinity 



2. Flow Inside the Circular Tube, Constant Heat Flux Case 

The graph of the nusselt number vs. x for the circular tube, constant heat flux flow 



was plotted for various Peclet numbers in figure 13 where x"^ is = 



2x1 d 
Re^Pr 



CIRCULAR TUBE CONSTANT WALL TEMP. FLOW. NUSSELT NUMBER 




Figure 13. The nusselt number distribution for circular tube constant heat flux flow 



The nusselt number is independent from the Peclet number for this case, and approaches 
to 4.364 where x goes to infinity. 
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3. Flow Between the Parallel Plates, Constant Wall Temperature 



Case 



The graph of the nusselt number vs. x for the parallel plates, constant wall 
temperature flow was plotted for various Peclet numbers in figure 14 where is 



2x!D^ 

RePr 



, Dh is 4h. 



PARALLEL PLATES CONSTANT WALL TEMP. FLOW. NUSSELT NUMBER 




Figure 14. The nusselt number distribution for parallel plates constant wall temperature flow 



The nusselt numbers approach to the values in those of table 2 when x gets larger as seen 
in figure 12. 



Pe 


0.5 


1 


2 


3 


5 


10 


50 100 


Nu. no 


8.0588 


8.0059 


7.9168 


7.8464 


7.7471 


7.6303 


7.5457 7.542 


Table 2. T 


le nusselt number for various Peclet numbers for parallel plates constant wall temperature case 



as X goes to infinity 
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^ Flow Between the Parallel Plates, Constant Heat Flux Case 

The graph of the nusselt number vs. x for the parallel plates, constant heat flux 



+ .+ 2x/Z);, 



flow was plotted for various Peclet numbers in figure 15 where x"^ is = 



RePr 



,Du is 



4h. 



PARALLEL PLATES CONSTANT HEAT FLUX FLOW, NUSSELT NUMBER 




Figure 15. The nusselt number distribution for parallel plates constant heat flux flow 



The nusselt number is independent from the Peclet number for this case, and 
approaches to 8.235 where x goes to infinity. 

C. HEAT ABSORPTION UPTO THE ENTRANCE 

The energy change from the entrance to x=0 where the heating starts has been 
investigated and plotted for fully developed laminar flow inside the circular tube and parallel 
plates for constant wall temperature and heat flux cases in figures 16-19. The Peclet numbers 
0.5, 1, 2, 3, 5, 10, 50, and 100 were used for the plots to see how the Peclet number, Pe^RePr, 
affects the heat absorption of the fluid. As it is seen from the graphs, the heat absorption 
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logarithmically increases as the Peclet number decreases. Therefore axial conduction effects or 
heat absorption upto the entrance where the heating starts has to be considered for low Peclet 
numbers. 

1. Flow Inside the Circular Tube, Constant Wall Temperature 

Case 

The graph of the heat absorption, E* vs. Peclet number for the circular tube, 
constant wall temperature flow was plotted in figure 16 where the dimensionless 
absorbed energy E* is £* = _ ^absorbed 

{to -tc)PCprl 




Figure 16. The energy absoption upto the entrance for circular tube constant wall temperature flow 



2. Flow Inside the Circular Tube, Constant Heat Flux Case 

The graph of the heat absorption, E* vs. Peclet number for the circular tube, 
constant heat flux flow was plotted in figure 1 7 where the dimensionless absorbed energy 



E‘ is £‘ = — 



' absorbed 



q'oPCpfolk' 
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Figure 17. The energy absoption upto the entrance for circular tube constant heat flux flow 



3. Flow Between the Parallel Plates, Constant Wall Temperature 



Case 



The graph of the heat absorption, E vs. Peclet number for the parallel plates, constant 
wall temperature flow was plotted in figure 1 8 where the dimensionless absorbed energy E* is 



E = 



“^absorbed 







Figure 18. The energy absoption upto the entrance for parallel plates constant wall temperature flow 
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4. Flow Between the Parallel Plates, Constant Heat Flux Case 

The graph of the heat absorption, E* vs. Peclet number for the parallel plates, constant 



heat flux flow was plotted in figure 19 where the dimensionless absorbed energy E’ is 



£' = 



^ absorbed 

qlbh^/a‘ 




Figure 19. The energy absorption upto the entrance for parallel plates constant heat flux flow 



D. AXIAL CONDUCTION EFFECTS INCLUDING THE WALL 
CONDUCTION 

The axial conduction effects for circular tube constant wall temperature case has 
been investigated including the wall conduction effects. The Peclet numbers 1, 3, 5, 
outer/inner radius ratios of 1.2, 2 and fluid/wall thermal conductivity ratios of 1 and 1/10 
were used to see how these properties affect the temperature, nusselt number and heat 
absorption. The theta bulk mean temperature, nusselt number variations were plotted vs. 
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x”" and the heat absorption were quantified as a function of Peclet number, thermal 
conductivity and radius ratios. The results are seen on the following pages. 

1. Bulk Mean Temperature Variations 

The ^c/ vs. graphs are seen in figures 20. a and 20.b for different Peclet 

^e/.x=0 



numbers, thermal conductivity and wall thickness. Since the bulk mean temperature 
values are very small about x'"=0, the normalized bulk mean temperature on the centerline 
is a good way to see the change at the bulk temperature relative to the centerline 

+ ^ • /I V~^o) t + • + 2x1 d 

temperature at x =0, where the bulk temperature is ^ ^ and x is x = . 

(l.-O RePr 

The bulk mean temperature vs. x'" graphs for the Peclet numbers 1,3, and 5, and the 
thermal conductivity ratios of 1, 1/10 and the radius ratios of 1.2 and 2 are seen in figure 



21 where the bulk mean temperature is 0^ = 







Figure 20.a The normalized bulk mean temperature with and without wall 
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Figure 20.b The normalized energy absorption upto the entrance for parallel plates constant heat flux flow 

The theta bulk mean values with the wall conduction case is closer to the steady 
state values than the no wall conduction case as it is seen in figure 21. The bulk mean 
temperature go to zero as x goes to plus infinity. This means that the temperature of the 
fluid including the wall conduction increases faster than the temperature of the fluid 
neglecting the wall conduction. The main reason of the relative increase at the 
temperature including the wall conduction case is the heat that is transferred by the axial 
conduction of the wall opposite to the flow direction. The heat moves to the opposite of 
the flow direction and is absorbed by the fluid. So the fluid temperature increases much 
more than no wall case. 




Figure 21. The bulk mean temperature with and without wall conduction 
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2. Nusselt Number vs. X 



The Nusselt numbers vs. x+ graphs are seen in figure 22 for different Peclet numbers, 

2jc / d 

thermal conductivity and wall thickness where is x”" = . The Nusselt number 

RePr 

distribution is very different from the no wall conduction case as it is seen in the figure. 
As the wall gets thinner and the thermal conductivity increases, the values of the Nusselt 
number gets closer to the case of no wall case. h=l/10 and k=1.2 is close to the no wall 
conduction case, so the values of the Nusselt number for these curves are close to the 
values of no wall case. 




Figure 22. The nusselt number distribution including the wall conduction 
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3. Heat Absorption vs. X 

The graphs of the heat absorption, E* vs. Peclet numbers for the conductivity ratios of 
1,1/10 and radius ratios of 1.2 and 2 were plotted in figure 23 where the dimensionless absorbed 

energy E* is £* _ ^ahsoM . As it is seen from the graph. As the wall's conductivity or 

(/<, -tc)PCp>'! 

thickness increases, the heat absorption before x=0 significantly increases. The wall 
conduction case has to be added to the problem for the microheat exchanger problems. 
Because the relative thickness of the wall is very large as it is compared to the inner 
diameter of the tube. The values of h=1.2 or h=2 are realistic for the microheat machines, 
because the diameter of the tubes are very small for these heat exchanger. The heat 
absorption for kf/ks=l/10 and h=1.2 case is more than ten times of the case of no wall. If 
the wall thickness would be h=2, then the heat absorption would be much greater. 




Figure 23. The energy absorption upto the entrance with and without the axial wall conduction 
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VI. CONCLUSIONS 



Axial conduction effects in duct flows are important and have to be considered for 
cases of Pe<10. The theta bulk mean temperature and Nusselt number variations vs. 
downstream distance were plotted as seen in the results. It is obviously seen the increase 
in the temperature against the flow direction from x=0 as the Peclet number decreases. 
Heat absorption by the fluid in the insulated region was also quantified by using the non- 
dimensional parameter E*. The heat absorption increases logarithmically as the Peclet 
number decreases, as seen in the heat absorption graphs. When the wall conduction 
effects are included in the problem, it was calculated that most of the heat absorption by 
the fluid is because of the wall conductivity as the wall conductivity gets larger. The same 
effect is also seen with increasing wall thickness. The heat absorption by the fluid gets 
larger when the wall gets thicker. It is very important to include the wall conduction in 
the problem as the wall gets thicker and the thermal conductivity of the wall gets larger. 
With large wall conductivity and wall thickness, the fluid temperature increases 
significantly before the heating starts. The current interests in axial conduction effects are 
because of the advent of micro heat exchanger technology with microtubes where the 
Reynolds number is very small. 
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APPENDIX - A SAMPLE PROGRAM 



The Program for Circular Tube, Constant Wall Temperature Case 

This program solves the finite difference equations for the circular tube, constant 
wall temperature case by using the Gauss Seidel Iteration Method with the overrelaxation 
parameter w=1.8. The program solves the equations doing the following steps : 

1 . The matrix is formed by the program, the first row represents the region which 
touches the tube. The last row is the centerline. The first column is the condition, where x 
goes to minus infinity. The last column represents the boundary condition as x goes to 
plus infinity. 

2. The boundary conditions are applied to the problem. The centerline boundary 
condition, the boundary conditions at r=ro, the boundary conditions as x goes to minus 
and plus infinity are applied to the last row, first row, first column and the last column 
respectively. 

3. The iterations are made calculating the values at all the grid points. As it is seen 
in the program, the odd rows are calculated first, the even rows later. The convergence of 
the problem is better in this case. 

4. The difference between the error and the boundary conditions for the plus and 
minus infinity is checked. If the boundary conditions and the difference between the 
iterations are reasonable, the program is stopped and the data is processed. The max. 
difference between the iterations is less than lE-6. The max difference between the 
values of the horizontal grid points which represents minus and plus infinity is less than 
lE-6. 



71 



% This problem solves the problem for circular tube 
% constant wall temperature case 
% 

% Programmed by Ibrahim Girgin 
% 



format long 

clear 
start=l ; 

iteration_number=300000 ; 

% Nx is the column number of the matrix 

% Ny is the row number of the matrix 

% X zero is the column number which represents x=0 

Nx=2501 

Ny=101 

x_zero=1501 

% For Peclet number is 1 



Pe=l 



% For the relaxation parameter 1.8 



w=1.8 

save iteration=100 



% hy is the vertical distance between the grids 
% hx is the horizontal distance between the grids 

hy=l/(Ny-l); 
hx=hy ; 

% teta is the matrix which holds the grid points 
% initialization of the teta 



teta=ones (Ny,Nx) ; 



% Apply the boundary condition for r=Ro, x>=0 



teta ( : , x_zero :Nx) =zeros (Ny, length (x_zero :Nx) ) ; 



% e .g. 

% the initial 
% value of 
% theta 
% 

% 



j=Ny >1 1 

1 1 
1 1 
1 1 
1 1 
j=l >1 1 



1 

1 

1 

1 

1 

1 



1 

1 

1 

1 

1 

1 



1 

1 

1 

1 

1 

1 



1 

1 

1 

1 

1 

1 



10 0 
10 0 
10 0 
10 0 
10 0 
10 0 



0 0 0 0 0 
0 0 0 0 0 
0 0 0 0 0 
0 0 0 0 0 
0 0 0 0 0 
0 0 0 0 0 



i=l 



i=x zero 



The values for this theta matrix is : 



i=Nx 
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Nx=14 Ny=6 x_zero=8 hx=0.2 hy=0.2 

Pe, w, f ile_name, save iteration is what is to be given 



% The variables that are being used inside the loop 

varl=Pe^2*hx^2 / 
var2=2/hy^2 ; 
var3=l/2/hx; 
var4=l/hy"2; 
var5=l/varl-var3 ; 
var6=l/varl+var3 ; 
var 7=1/2 /hy; 
var8=2/varl+4*var4 ; 
var9=2/varl+var2 ; 
varl0=2*hx; 

% The iteration starts here. For every value of kkk, one 
% iteration is made 

for kkk=start : iteration number 



for max err=teta; 



for i=2:Nx-l 
first 



for j=l:2:Ny 
y= (j -1) *hy; 
if j==l 



^ do the iteration for all the columns except the 
% and the last columns 

% do the iteration for all the odd rows 
% calculate the y, 0<= y <=l 

% do the iteration for radius=0 (for the center) 



teta (Ny, i) = (teta (Ny, i+1) * (var 5) +teta (Ny, i-l) * (var 6) + ... 
4*teta (Ny-1, i) *var4) /var8 ; 



elseif j==Ny 



do the iteration for radius=R 



if i<x zero 



if i<x zero & for radius R, it is insulated 



teta (1 , i) = ( (teta (1 , i+1) +teta (1 , i- 
1) ) / (varl) +teta (2 , i) *var2) /var9 ; 



end 



else 
and r=R 



% do the iteration for all the odd rows except r=0 



teta (Ny- j +1 , i) =teta (Ny- j +1 , i) +w* ( (teta (Ny- j +1 , i+1) * (var5 + 
y^2/varl0) . . . 



73 



+teta {Ny-j+1, i-1) * (var6 
var7/y) + . . . 

teta (Ny- j +2 , i ) * ( var4 - 



- y^2/varl0) +teta (Ny- j , i) * (var4 
var7/y) ) /var9-teta (Ny-j+1, i) ) ; 



+ 



end 



end 

for j=2:2: (Ny-1) % do the iteration for all the even rows 

y= ( j -1) *hy; 

teta (Ny-j+1, i) =teta (Ny-j+1, i) +w* ( (teta (Ny-j+1, i+1) * (var5 + 
y^2/varl0) . . . 

+teta (Ny-j+1, i-1) * (var6 - y^2/varl0) +teta (Ny- j , i) * (var4 + var7/y)+ 
teta (Ny- j+2, i) * (var4 - var7/y) ) /var9-teta (Ny-j+1, i) ) ; 



end 



end % End of the iteration 



Save the variables at every "save_iteration” 
e.g. if the save_iteration equals to 20, then save 

the variables at when kkk=20 , 40 , 60 , 80 , . . . so on. 



if fix (kkk/save_iteration) ==kkk/save_iteration; 

max_error= (max(max(abs ( teta-for_max_err) ) ) ) 
clear for_max_err 

save forpel teta Nx Ny x_zero w Pe hx hy max_error kkk 
save_iteration start y 



end 

end 



74 



LIST OF REFERENCES 



[1] Frank P. Incropera, David P. Dewitt, Introduction to Heat Transfer, John Wiley & 
Sons, 1996, Canada. 

[2] W. M. Kays, M. E. Crawford, Convective Heat and Mass Transfer, McGraw-Hill, Inc. 
1993. 

[3] Adrian Bejan, Convection Heat Transfer, John Wiley & Sons, 1996, Canada. 

[4] L.Graetz, Uber die Warmeleitungsfahigkeit von Flussigkeiten (On the thermal 
Conductivity of liquids) Part 1. Ann. Phys. Chem., Vol. 18, pp.79-94, 1883,Part2; 
Ann.Phys. Chem., Vol. 25, pp. 337-357,1885. 

[5] M.L.Michelsen and John Villadsen, “The Graetz Problem With Axial Heat 
Conduction,” InternationalJournal of Heat Mass Transfer, Vol.l7, pp. 1391-1402, 1974. 

[6] X. Yin and H. H. Bau, “The Conjugate Graetz Problem With Axial Conduction,” 
Journal of Heat Transfer, Vol. 11 8, pp.482-485. May 1996. 

[7] M. L. James, G. M. Smith, J. C. Wolford, Applied Numerical Methods for Digital 
Computation, Fourth Edition, HarperColIinsCollegePublishers, 1993. 

[8] George E. Forsythe and Wolfgang R. Wasov, Finite Difference Methods for Partial 
Differential Equations, Wiley, pp.256, 1960, New York. 



75 



76 



INITIAL DISTRIBUTION LIST 



No. Copies 

1 . Defense Technical Information Center 2 

8725 John J. Kingman Rd., STE 0944 

Ft. Belvoir, VA 22060-6218 

2. Dudley Knox Library 2 

Naval Postgraduate School 

41 1 Dyer Rd. 

Monterey, CA 93943-5101 

3. Chairman, Code ME 1 

Department of Mechanical Engineering 

Naval Postgraduate School 
Monterey, CA 93943-5000 

4. Professor Ashok Gopinath, Code ME/GK 2 

Department of Mechanical Engineering 

Naval Postgraduate School 
Monterey, CA 93943-5000 

5. Naval/Mechanical Engineering Curricular OfBce, Code 34 1 

Naval Postgraduate School 

Monterey, CA 93943-5000 

6. Deniz Kuwetleri Komutanligi 1 

Personel Daire Baskanligi 

Bakanliklar, Ankara, TURKEY 06100 

7. Deniz Harp Okulu Komutanigi Kutuphanesi 2 

Tuzla, Istanbul, TURKEY 81704 

8. Ibrahim Girgin 3 

Alti Eylul Mah. Aksu Sok. No;25 

45450 Ahmetli, Manisa, TURKEY 

9. Istanbul Teknik Universitesi Makina Muh.Bl. 1 

Ratip Berber Kutuphanesi 

Maslak, Istanbul, TURKEY 



77 



15 

10/99 



4B3NPG 
TH 

22527-200 



r 





4 4 ,. 



• 'V 










r 






L.. 






> 







' A 
** 



W' 



>■ 




















'3 



m 







